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Abstract 

We demonstrate that it is possible to use the balanced homodyning with 
array detectors to measure the quantum state of correlated two-mode sig- 
nal field. We show the applicability of the method to fields with complex 
mode functions, thus generalizing the work of Beck(Phys. Rev. Lett. 84, 
5748 (2000)) in several important ways. We further establish that, under 
suitable conditions, array detector measurements from one of the two ouputs 
is sufficient to determine the quantum state of signal. We show the power 
of the method by reconstructing a truncated Perelomov state which exhibits 

complicated structure in the joint probability density for the quadratures. 
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I. INTRODUCTION 



The subject of quantum state tomography (QST) has been of great interest in the recent 
years after the experimental reconstruction of the complete wave function of the vacuum 
and squeezed states by balanced homodyne (BHD) method Much progress has been 

made by way of theoretical exploration of alternate schemes of QST like the self-homodyne 
tomography and their experimental verification [Q. In essence, the BHD method combines 
a strong local oscillator (LO) with the signal in a 50:50 beam-splitter (BS). The two outputs 
of the BS are measured by two photodetectors and electronically subtracted. The resultant 
quantity is directly proportional to the rotated quadrature [aexp(— i<p) + a) exp(i0)]/\/2 of 
the signal field and the angle is fixed by the LO. By performing the experiment many 
times for a given <p and repeating it for various values of <fi in the range of < < 7r, the 
probability density p(x, <j)) for the quadrature rotated through <p can be measured. From 
the measured probability density the Wigner function for the signal can be constructed 
and hence the elements of the density matrix [||||. The serious problems that may arise 
in the numerical reconstruction of density matrix elements and the methods to avoid such 
pitfalls are discussed in the review article by Welsch et al |J and in Ref. |7[]. The BHD 
can be extended p|-p^] to multimode signals and to measure other distributions like the 
positive P-distribution. The BHD is for optical fields and methods have been developed for 
reconstructing the vibrational state of molecules |13[] and spin systems |14|]- For measruing 



the quantum states of the modes of a cavity, atoms with properly chosen energy levels can 
be used as probes [|l^ . 

A serious drawback of the BHD using single detectors is the requirement for mode match- 
ing the signal and LO mode functions. A quantitative measure of mode matching is the 
overlap between the signal and LO mode functions and the efficiency of the BHD scheme is 
directly proportional to it, 

efficiency oc / U* gnal (x)Ui J o(x)dx. (1) 
Here U S i gna \ and ?7lo are the respective mode functions of the signal and LO. If there is 



perfect mode matching, i.e., the overlap integral is unity, the BHD is very efficient. In the 
case of signal and the LO modes being orthogonal to each other, the BHD simply fails to 
determine the quantum state of the signal. 

An ingenious way of circumventing the problem of mode matching was suggested by 
Beck [TEfl. The key idea is to replace the single detectors at the BS output ports by array 



detectors. Beck considers array detectors with pixels small enough such that the mode 
functions of the LO and signal over any pixel is constant. The output from a given pixel 
labeled (j,f) of one of the array detectors is given by 

Ni{j,j'A) = 5j ^[^ + J2al l a m U:(j,j')U m (jJ') 

H — / = Y^[ & mexp(-i<l))U m (j,j') +H.c]\. (2) 

D-xDy m 



Here the signal is taken to be multimode and the LO is assumed to be a single mode coherent 
state |/3exp i?i ) (/3 real). The summation is over the various modes of the signal field. The 
function U m (j,j') is the value of the rath mode function of the signal over the pixel (j, j'). 
The operators and a m are the creation and annihilation operators for the rath mode of 
the signal. The mode function for the LO is taken to be the constant , 1 over the entire 

y/D x D y 

detector surface. The dimension of each pixel is Sx x 8y and that of the detector is D x x D y 
(=N5x x N8y). The constant c is the speed of light vacuum, T is the duration of counting 
and L is the longitudinal quantization length. Similarly, the output of the corresponding 
pixel at the other output port of BS is given by 



ALi L) x L)y nm 

Y^[a m exp(-i(p)U m (j,j') + H.c]\. (3) 



D X Dy 
D X Dy -III 

On taking the difference of Ni(j,f, 0) and N 2 (j,j', 0), it is seen that 

N d {j,j'A) = -4^=PT,^mexp(- i( p)U m (j,f)+H.c}, (4) 

D x D y m 

where we have set L = cT. If one of the mode functions, say Ui, is real, then the above 



expression yields, on summing over all the pixel indices and using orthonormal relationship 
among the modes, 

^N d (j,j',(p)Ui(j,j') = 1 /3[aiexp(-i<f>) + of exp(z </>)]. (5) 



It is to be noted that the mode matching factor does not enter in the expression relating 
the measured photocurrent difference and the rotated quadrature. This is of great practical 
utility as mode matching is extremely difficult in an actual experiment. The above expression 
for quadrature works for any mode function that is real and hence the quantum state of all 
such modes can be determined. However, their joint state cannot be estimated by the 
method as proposed by Beck. But correlated quantum states are of prime importance for 
conceptual understanding of quantum theory as well as applications |17|.lq|. Hence, it is all 
the more essential that the array detector scheme be extended to correlated multimode fields 
whose mode functions need not be real. 

In the present work we extend the Beck's method to two-mode correlated states by 
interposing a linear device which introduces a variable relative phase between the two modes. 
The stringent condition to have real mode functions is relaxed by making measurements 
with the LO phase varying from zero to 3| and combining the measurements to get the 
value of the quadratures roatated through zero to tt as required for state reconstruction. 
Further, we establish that one of the outputs of BS as measured by an array detector is 
sufficient to measure the quantum state. This paper is organized as follows. In Section II 
the extension of Beck's work to two-mode correlated state is given. Section III of the paper 
contains the results of Monte Carlo simulations to demonstrate the power of the method 
to construct joint probability densities with complicated structures. In Section IV of the 
paper we explicitly show that for mode functions satisfying certain conditions, a single array 
detector is sufficient to measure the quantum state. We have provided an appendix to 
make transparent the arguments given Section II for realxing the requirement for real mode 
functions. 
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II. JOINT QUANTUM STATE OF TWO-MODE FIELD 



In the present section a method of estimating the quantum state of a correlated two 
mode signal is given. The basic requirement is to construct the joint probability density 
P(x\, <j>i, X2, 4>2) f° r the two quadratures rotated through 0i and 2 respectively. It is, 
therefore, essential to rotate the two quadratures independently. If the signal modes are 
spatially seperable, they can be mixed with two independent LOs (usually made from a 
single source by using a BS) and measurements are carried out using two BHD arrangements 
M. In the other case, the signal is passed through a medium which mixes the modes with 
a relative phase. The mixed modes are taken to be spatially seperable [g] and one of the 
seperated modes is used in the input port of BS. This requires only one LO and one BHD 
arrangment. In both the cases, however, spatial seperation of modes is required at some 
stage and mode matching is essential. In the present section we discuss how to use the BHD 
with array detectors in the latter method. The mixing of modes is achieved by using a linear 
device which introduces a relative phase between the two modes. For instance, the linear 
device could be a medium in which the two modes interact via the Hamiltonian given by 

H int oc djd 2 exp(iO) + H.c (6) 

The field operators (a' v a! 2 ) at the output of the linear device can be written as 

a'i = cosz/di — isin z/exp(— i6) d 2 , (7) 
a' 2 = cosz/d 2 — isin uexp(i9) di, (8) 

where v is a constant determined by the length of the linear device, strength of interaction, 
etc. The device mixes the two modes of the signal. For later use we define the rotated 
quadratures for the transformed field: 

X-(0,^,,) = ^ eXp( - #) + a " eXpW , jfe = l,2 (9) 

V2 

The positive frequency part of the electric field operator entering the signal port of BS is 
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EM'(x,y) 



2-kTioj 



1/2 



[U x (x,y)a' x + U 2 (x,y)a 2 ) 



(10) 



The electric field operator depends on 0,6* and v in addition to the spatial coordinates (x, y). 
However, in the sequel we include the dependence explicitly in the expressions for difference 
counts and other derived quantities. With the field operators given by Eq.|7] entering the 
signal port of BS and the LO phase fixed at 0, the difference count N d (j,f, 0, 9, v) is 

5x5y 



N d {j,j',cl>,e,v)=(3 



D X Dy 



exp(-20)[£/i(j, j')a[ + U 2 (j,j)a 2 ] 
+ exp(#) [U* (j, 3 ')a ] ( + U* 2 (j, f)at\ 



(11) 



and the difference count with the LO phase rotated further by it/ 2 is 



N d ( 3 ,f,<f ) + ^e,u) = i (3 () ' r,)! ' 



-exp(-^0)[t/ 1 (j,j')a / 1 + U 2 (j,f)af 2 ) 



(12) 



The quantities N d (j, j', 0, 9, v) and N d (j, f, <f> + %,9,is) are measured in the experiment. 
From these measured quantities we construct 



R(j, f, 0, 0, „) = i JM>L{N d (j, f, 0, 9, v) 



2(3 V 2 
-iJ\ 
5x5y 



71 



iN d (j,j',(t ) + -,9,u)}, 



exp(i<j))[a^U*(j,f) + a%XJ* 2 (j, j')}. 



V2 



(13) 
(14) 



Multiplying R(J, j', 0, 9, v) by Uk{j,j') (k = 1,2) and summing over (j,f) yields 



Y, R {3,3'A,0,v)Uk{j,j') = ^=exp( 
j,f v2 
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A- ' 



(15) 



The quadrature X' fc (0) can now be written in terms R(j,f, 0, 0, v) and its conjugate as 



x'M = Y\R{3,?A,e,v)Uk{3,j') + c.c]. 



(16) 



The equations determine the mixed quadratures in terms of measured difference counts. 
Using Eqs. J7[||, the quadratures for the mixed mode can be written in terms of the signal 
quadratures as 



X\((f>,6,v) = cosvX 1 (4>)+smvX 2 (4> + 6 + ^), (17) 
X' 2 ((f>,6,v) = cosz/X 2 (0) + sinz/X 1 (0-# + ^). (18) 

If measurements are carried out for two values of u, say V\ and u 2 , then the signal quadratures 
can be evaluated from Eqs.(|iT|)-(|i"8|) and the resulting expressions are: 

x 2 (4>) 
li(0-0 + |) 

71 
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V 2 - 
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>i - 
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^2) 


sin( 


V 2 - 
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^1) 



sin V\X' 2 {Q, 


0,^2) 


- sinz/ 2 X' 2 (#, 


0,^i)], 


(19) 


COS UiX'z(0, 


0,^2) 


— cos z/ 2 X' 2 (#. 


0,^i)], 


(20) 


sin 


0,^2) 


— sin v 2 X \(9, 


0,^i)], 


(21) 


cos 


0,^2) 


— cos v 2 X\{9, 


0,i/i)]. 


(22) 



^2(0+^ . . , 

2 sin(^ — z/ 2j ! 

The first two expressions yield the values of the quadratures X 2 and X\ rotated through 
angles <ft and <fi — 9 + | respectively. The last two yield the value of two quadratures rotated 
through <pi = <p and 2 = + 9 + |. By changing the LO phase and the phase # in the 
interaction Hamiltonian, the quadratures can be measured to construct the joint probability 
distribution P(x\, 0i, x 2 , 2 ). Note that there is no need to assume that the mode functions 
are real. Of course, the difference count measurements are to be carried out over a range 
of — Ztt/2 for the LO phase; in the case of real mode functions it is sufficient to measure 
the difference count when the LO phase is varied from — tt. From the measured joint 
probability density, the Wigner function for the two-mode state can be constructed which, 
in turn, can be used to determine the density matrix. The measured joint probability can 
also be used to construct the elements of the density matrix in two-mode number state basis 

The scheme suggested here is similar to what is described as the method of generalized 
rotation in phase space 0. However, in our scheme there is no need to spatially separate 
the signal modes at the output of the linear device. This is possible as mode matching is not 
essential in the present scheme. In fact, not separating them is useful to reduce the number 
of experimental runs. In one run of the experiment, the pairs {Xi(0),X 2 (0 + 9 + |)} and 
{Xi(<p — 9 + |),X 2 (0)} are simultaneously estimated. 
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III. MONTE CARLO SIMULATIONS 



In this section we study the experimental feasibility of the above method by Monte 
Carlo (MC) simulations with the signal in the two-mode Perelomov state P0j. This state is 
generated from the two-mode vacuum |0, 0) as follows: 



|C) = exp(Ca t 1 a 2 - C*aia 2 )|0,0). (23) 



The number state expansion for the above state is given by 



J oo 

10 = — t— ^[exp(-z 7 )tanhr]"in,n), (24) 
cosn r n=0 

and the joint probability density is 

p(xi,<p 1 ,X2,<h) = —£sexp[ j - J. (25) 

The constants r, A and B are related to the squeeze parameter (: 



C = rexp(-ij), 



A 
B 



|1 + tanhr exp[-i(0i + 2 + 7)]| 



1 — tanh 2 r 



|1 - tanhr exp[-z(0i + 2 + 
1 — tanh 2 r 



2 



2 



Random numbers distributed according to the distribution given in Eq. |25] were gener- 



ated from Gaussian distributed random numbers by von Neumann's rejection method [21 
These numbers were then used to generate the output of an actual experiment. In Fig. 1 
the result of the simulation to reconstruct the joint probability density is given for <fii = 7r/4, 
02 = tt/2, 7 = 7r/4 and r = 1.0 along with its contour plot . For comparison, the theoretical 
distribution is also given. We find that to reconstruct the probability density with reasonable 
accuracy, as many as 1.6 x 10 5 experimental runs would be required for each set of values of 
<pi and 02- I n Fig- 2 we present the simulated results for identical values for the parameters 
except 02 which is set equal to — 7r/4. In Fig. 1 the distribution is narrower compared to 
what is shown in Fig, 2. As expected, from the figures it is clear that the reconstruction is 
better, for a given number of experiments, if the probability density is narrower. 
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In order to see whether an experiment can capture the complicated spatial structures of 
probability density, the MC simulations were done for the signal in truncated Perelomov 
state. These are defined as 



| Cl ,c 2 ) = ci|0,0) + c 2 exp(-^)|l,l). (26) 

The superposition coefficients c\ and c 2 are taken to be real and they satisfy the normaliza- 
tion condition c\ + c| = 1. The joint probability density for the quadratures is 

p(x!, 4>i,x 2 , 4> 2 ) = [c\ + A<? 2 x 2 x 2 + 4:X 1 x 2 cic 2 cos(0i + 2 + 5)]. (27) 

7T 

The result of the MC simulation for the truncated Perelomov state is depicted in Fig. 3. 
We have used c\ = c 2 = l/\/2, 5 = 7r/8, <f>i = ir/4 and <p 2 = n/A. It requires 2 x 10 5 points 
to generate the pattern shown. 



IV. QST USING A SINGLE ARRAY DETECTOR 

In this section we describe how to measure the joint probability density of a two-mode 
state using an array detector in one of the output ports of BS. No measurements are made 
at the other output. In the discussions to follow the pixels are labeled by a single index 
instead of two indices. We first introduce some notations to simplify the expressions: 

Ui(k) : value of the mode function U\ on the pixel labeled k, 
U 2 {k) : value of the mode function U 2 on the pixel labeled k, 
n d (k) : difference count between the pixels labeled 1 and k (k — 2, 3. ..9) 
V T : transpose of vector V whose elements are 

{aj'%, a 2 a' 2 , a|'a 2 , o! x o) 2 \ exp(— ifya^, exp(i0)a{', 
exp(— i(f))a' 2 , exp(i0)a 2 '} 
M : an 8 x 8 matrix whose (k — l)th (k = 2, ..9) row is given by 
{£/?(l) -!/?(*), Ul{l)-U 2 2 {k), C/*(l)[/ 2 (l) - U*(k)U 2 (k), 
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*7i(l)£/ 2 *(l) - U^QUZik), Sp^l) - U x (k)}, 
S[U*(l)-U;(k)}, S[U 2 (1) - U 2 (k)}, s[u;(i)-ut(k)]h 

p 



where S 



D x D y 



The expression for rid{k) (k — 2, 3 • • • 9) is 

nd{k) = Counts in pixel 1 - Counts in pixel k 

5x5y I 



2 



[ E a]:a' m [U* n (l)U m (l) - K(k)U m (k)) 



n,m=l 



+ -J= £[exp(-i0)a' m [[/ m (l) - £/ ro (A:)] + H.c)] (28) 



D-r.Dii ) 



= ^EW (29) 
z i=i 

The elements of the matrix are various combinations of the mode functions U\ and U 2 over 
a selected set of nine pixels of the array detector. If the elements of M are such that the 
matrix is invertible, the elements of the vector V can be determined from the equation 

9 

V j = Y / (M- 1 ) j>k . 1 n d (k) j = l,2---8. (30) 

k=2 

Note that the vector V has as its elements the creation and annihilation operators of the two 
modes and their quadratic forms. The above equation then implies that the elements of V 
are determined from the measured quantities nd(k). It is evident that the mixed quadratures 
can be determined as 

X\{4,0,u) = yi±p., (31) 
X' 2 fae,v) = ?l±p-. (32) 

Once the mixed quadratures are estimated, the signal mode qudratures are obtained using 
Eqs. |T7| - |I8| . It is clear that if we can choose nine pixels so that the matrix M is invertible, 
then the joint probability density of two-mode states can be determined. Despite the fact 
that only one of the output ports of BS is used, the LO fluctuations are eliminated as in the 
case of BHD. Another interesting feature is that there is no need to assume that the mode 
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functions are real and the LO phase need not be changed beyond the usual range — n. 
The method can be extended to more than two modes. But the requirements on the mode 
functions will become more stringent. 



V. SUMMARY 

The use of array detectors in BHD eliminates the need for difficult-to-achieve mode- 
matching. For two-mode signals, a linear device to introduce a relative phase between two 
modes can be used in the signal port of the BS. By varying the LO phase from zero to ^ and 
the relative phase between the modes, the joint probability density for the quadratures of the 
two modes can be experimentally determined. The QST measurements can carried out with 
one array detector instead of detectors at both the output ports of BS. The method retains 
the advantage of BHD in the sense that the LO fluctuations do not affect the measurement. 
However, the success of the method depends of the shape of the mode functions. 



APPENDIX A: QST OF A SINGLE, COMPLEX MODE SIGNAL 

In this appendix we show how to measure the quantum state of a single, complex mode 
signal. This appendix is provided to make the arguments of Section || more transparent. 
In Beck's method it is required to have the mode function real for determining its quantum 
state. If the mode function is not real, then it cannot be factored out in the RHS of Eq. (|J) 
and hence the measured difference count cannot be related to the quadrature. Note that 
this problem does not arise in the case of the conventional BHD using single detectors. We 
set Nd(j,j',4>) as Nd(<f>) and U(j,f) as U to avoid lengthy expressions. Hence, Eq. (|) is 
rewritten as 

N d (cj)) = 5x5y (3[U exp(-i(f))a + U* exp(i(f))a ] ]. (Al) 
y D x D y 

With the LO phase rotated through 0+5, the difference count is 
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N d ((j)+ = i y p[tf exp(i<j>)U* - U exp(-i<p)a]. (A2) 

2 ^ D X Dy 

Comibining the expressions for Nd(4>) and Nd(<f> + § ) we get 

W) - tN d (<f) +\) = U* exp(z0)a t (A3) 

and 

iV d (0) + iJV d (0 + |) = Uexp{-i<f>)a. (A4) 

Having measured the operators exp(— i0)a and its adjoint, the quadrature can be estimated 
at once as it is a linear combination of the two operators. To achieve this we multiply the last 
two equations by U and U* respectively and use the normalization relation 5x5y U*U = 
1 to yield 



x{4>) = EWW + ir) + iNM + |)(tr - [/)]. (as) 

As asserted in the beginning of the appendix, there is no need to assume the mode functions 
to be real in the above expression. When U is real the second term on the RHS vanishes 
and we recover, as expected, Eq.^|. 
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FIGURES 

Thoeretical 



Simulated 




FIG. 1. Joint probability density for two-mode Perelomov state. The parameter values are 
(p\ = 7r/4, (j>2 = vr/2, 7 = 7r/4 and r = 1.0. 
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Thoeretical 



Simulated 




FIG. 2. Joint probability density for two- mode Perelomov state with with the parameters taking 
the same values as in Fig. 1 except 4>2 which is set equal to — 7r/4. 
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Thoeretical 



Simulated 




FIG. 3. Theoretical and simulated joint probability densities for the truncated Perelomov 
state. The parameter values are c\ = C2 = l/y/2 and 5 = tt/8 
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